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ABSTRACT 


In  this  investigation,  a  method  for  the  finite  rotation  and  large  deformation  analysis  of  plates  is 
presented.  The  method,  which  is  based  on  the  absolute  nodal  coordinate  formulation,  leads  to  an 
isoparametric  plate  element  capable  of  representing  exact  rigid  body  motion.  In  this  method, 
continuity  conditions  on  all  the  displacement  gradients  are  imposed.  Therefore,  non-smoothness  of 
the  plate  mid-surface  at  the  nodal  points  is  avoided.  Unlike  other  existing  finite  element 
formulations  that  lead  to  a  highly  nonlinear  inertial  forces  for  three-dimensional  elements,  the 
proposed  formulation  leads  to  a  constant  mass  matrix,  and  as  a  result,  the  centrifugal  and  Coriolis 
forces  are  identically  equal  to  zero.  Furthermore,  the  method  relaxes  some  of  the  assumptions  used 
in  the  classical  and  Mindlin  plate  theories.  By  using  a  general  continuum  mechanics  approach,  a 
relatively  simple  expression  for  the  elastic  forces  is  obtained.  By  developing  such  an  isoparametric 
plate  element,  the  development  of  three-dimensional  shell  elements  becomes  straightforward. 
Numerical  results  are  presented  in  order  to  demonstrate  the  use  of  the  proposed  method  in  the  large 
rotating  and  deformation  analysis  of  plates. 
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1. 


INTRODUCTION 


Finite  element  formulations,  including  incremental  methods  [1,  2,  10]  and  large  rotation  vector 
formulations  [15,  16]  are  widely  and  successfully  used  in  the  analysis  of  structural  systems.  The 
applications  of  these  formulations  to  the  analysis  of  large  deformation  of  multibody  systems  have 
been  a  subject  of  interest  in  many  investigations.  The  incremental  methods  lead  to  a  linearization  of 
the  rigid  body  kinematic  equations  when  conventional  finite  elements  are  used.  In  order  to 
overcome  this  problem,  several  large  rotation  vector  formulations  were  recently  proposed.  In  these 
formulations,  finite  rotation  parameters  are  used  as  nodal  coordinates.  Continuity  conditions  are 
imposed  on  the  nodal  displacements  and  the  rotation  parameters  at  the  element  interfaces.  The 
concerns  with  regard  to  the  use  of  the  large  rotation  vector  formulations  in  the  analysis  of  flexible 
multibody  systems  can  be  summarized  as  follows: 

1 .  Imposing  continuity  on  the  finite  rotations  at  the  nodal  points  does  not  guarantee  the  continuity 
of  the  displacement  gradients  at  these  points.  As  a  result,  the  centerline  or  the  mid-surface  of  the 
element  is  not  smooth.  The  obtained  solution  eventually  leads  to  errors  in  the  calculations  of  the 
elastic  forces  and  stresses  at  the  nodal  points.  This  is  one  of  the  reasons  that  most  solutions  for 
the  large  rotation  vector  formulations  are  based  on  incremental  procedures  and  require  an 
elaborate  updating  scheme  for  the  finite  rotation  parameters. 

2.  The  interpolation  of  finite  rotations  should  be  carefully  handled,  particularly  in  three- 
dimensional  applications.  For  example,  linear  interpolation  for  Euler  parameters  does  not  imply 
linear  interpolation  for  Euler  angles  [6]  or  Rodrigues  parameters.  Different  sets  of  orientation 
parameters  are  related  by  highly  nonlinear  equations  [13].  This  fact  explains  the  need  to  use  an 
accurate  updating  scheme  for  the  finite  rotation  parameters  in  large  rotation  vector  formulations 
in  order  to  obtain  reasonable  solutions.  In  fact,  in  many  of  the  currently  used  updating  schemes. 
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an  incremental  procedure  is  used  by  which  the  finite  rotation  parameters  are  incrementally 
updated  to  allow  exploiting  the  features  of  infinitesimal  rotations. 

3.  Most  finite  element  large  rotation  vector  formulations  cannot  be  used  to  develop  a  more  general 
beam  and  plate  theory.  In  the  excellent  work  of  Simo  [15],  a  coordinate  system  is  attached  to  the 
element  cross  section.  The  equations  of  motion  of  the  cross  section,  treated  as  rigid,  are 
developed  in  terms  of  a  set  of  generalized  coordinates.  These  generalized  coordinates  that 
include  finite  rotation  parameters  and  depend  on  the  spatial  coordinates  of  the  elements  are  then 
interpolated  using  polynomials.  By  using  this  procedure,  it  becomes  difficult  to  relax  the 
assumption  of  the  rigidity  of  the  cross  section  of  the  element.  This  assumption  is  used  in  the 
classical  beam  and  plate  theories  as  well  as  many  of  Timoshenko  beam  [5]  and  Mindlin  plate  [8] 
models. 

4.  The  finite  element  large  rotation  vector  formulations  (as  a  consequence  of  3)  do  not  lead  to  a 
constant  mass  matrix  in  the  case  of  three-dimensional  applications,  and  as  a  result,  the 
centrifugal  and  Coriolis  forces  are  not  equal  to  zero. 

Most  of  these  fundamental  concerns  can  be  adequately  addresses  using  the  absolute  nodal 
coordinate  formulation.  It  is  a  simple,  non- incremental  procedure  that  leads  to  a  constant  mass 
matrix  in  two-  and  three-dimensional  analysis.  It  leads  to  smooth  centerline  and  mid-surface  of  the 
elements  and  ensures  the  continuity  of  the  rotations  of  the  cross  section  as  well  as  the  displacement 
gradients  at  the  nodal  points.  It  can  also  be  used  to  systematically  relax  some  of  the  assumptions 
used  in  Euler-Bernoulli,  Timoshenko  beam  and  Mindlin  plate  models.  Furthermore,  the  absolute 
nodal  coordinate  formulation  does  not  require  the  interpolation  of  finite  rotations  and  leads  to  exact 
modeling  of  the  rigid  body  dynamics. 
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One  of  the  features  of  the  non-incremental  absolute  nodal  coordinate  formulation  is  that  a 


global  coordinate  system  is  used  for  the  definition  of  the  nodal  coordinates.  In  this  case  no 
coordinate  transformation  is  required.  The  nodal  coordinates  and  the  element  global  shape  function 
are  used  to  define  the  location  and  deformation  of  material  points  on  the  finite  element.  It  is, 
therefore,  required  that  the  global  shape  function  must  be  able  to  describe  the  rigid  body  motion  as 
well  as  the  deformation  of  the  finite  element.  Another  important  feature  of  the  formulation  is  that 
global  slopes  instead  of  angles  are  used  to  define  the  configuration  of  the  finite  element.  The  use  of 
the  infinitesimal  rotations  in  the  classical  finite  element  formulation,  as  previously  pointed  out, 
leads  to  linearization  of  the  rigid  body  kinematic  equations,  and  as  a  consequence  the  description  of 
rigid  body  displacements  may  not  be  exact  [12].  When  global  slopes  are  employed  in  the  absolute 
nodal  coordinate  formulation,  beam  and  plate  elements  become  isoparametric  while  in  the  classical 
finite  element  formulation  such  elements  are  considered  non-isoparametric. 

The  objective  of  this  study  is  to  develop  a  new  finite  plate  element  based  on  the  absolute 
nodal  coordinate  formulation.  In  this  investigation,  a  more  general  plate  element  is  developed  by 
relaxing  some  of  the  assumptions  of  the  classical  and  Mindlin  plate  theories  [17].  The  mass  matrix 
of  the  more  general  plate  element  remains  constant,  and  therefore,  the  centrifugal  and  Coriolis 
forces  are  equal  to  zero.  By  using  a  continuum  mechanics  approach  and  nonlinear  strain- 
displacement  relationships,  a  general  expression  for  the  elastic  forces  that  takes  into  account  the 
effect  of  shear  deformation  is  obtained.  Results  are  presented  in  this  study  in  order  to  demonstrate 
the  generality  of  the  new  plate  model. 

2.  DISPLACEMENT  FIELD 

The  choice  of  the  type  and  form  of  the  interpolation  functions  plays  an  important  role  in  developing 
a  new  type  of  finite  elements.  In  this  study,  interpolation  polynomial  functions  are  selected  because 


they  are  easy  to  integrate  and  differentiate.  The  more  complex  deformation  modes  within  one 
element  can  be  achieved  by  increasing  the  order  of  polynomials  that  define  the  displacement  field 
of  the  plate.  The  use  of  higher  order  elements  leads  to  a  substantial  decrease  in  the  number  of 
elements  required  to  develop  a  discretized  structure  model.  In  the  absolute  nodal  coordinate 
formulation,  the  displacement  field  must  include  the  rigid  body  modes  as  well  as  the  plate 
deformation  modes.  The  following  polynomials  are  assumed  to  define  the  displacement  field  of  the 
plate  element  used  in  this  investigation: 


r  = 


V 

r2 

/3_ 

2  2  2  2 
ciq  +  axx  +  a2y  +  a3z  +  a4xy  +  a5xz  +  a6yz  +  a7x  +  a%y  +  a<4x  y  +  awxy 

2  2  2  2 

b0  +  b\X  +  b2y  +  b3z  +  b4xy  +  b5xz  +  b6yz  +  b7x~  +  b$y  +  b^x  y  +  b\$xy 

2  2  2  2 

c0  +  cxx  +  c2y  +  c2z  +  c4xy  +  c5xz  +  c6yz  +  c7x  +  c8y  +  c9x  y  +  cl0xy 


3  3  3  3 

+  at  {xyz  +  a{ 2x  +  auy  +aX4x  y  +  ax 5xy 

+  bx  xxyz  +  bl2xJ  +  bl3yJ  +  b^  y  +  bX5xy^ 
3  3  3  3 

+  q  \xyz  +  Ci 2x  +  Cj 3y  +  Cj 4x  y  +  cl5xy 


(1) 


where  r  is  the  global  position  vector  of  an  arbitrary  point  P  on  the  plate,  x ,  y  and  z  are  the  spatial 
coordinates  of  the  plate  element  coordinates  and  a„  b,  and  c,  are  polynomial  coefficients  that  depend 
on  time  in  the  case  of  dynamic  problems.  The  spatial  coordinate  z  is  assumed  to  be  in  a  direction 
perpendicular  to  the  mid-surface  of  the  plate  in  the  undeformed  state,  while  the  spatial  coordinates  x 
and  y  are  measured  in  the  mid-surface  plane.  The  displacement  field  in  the  Eq.  1  is  expressed  in  the 
global  coordinate  system  and  therefore  the  order  of  the  three  polynomials  are  the  same.  It  can  be 
demonstrated  using  Eq.  1  that  the  assumed  displacement  field  couples  the  rigid  body  motion  and  the 
plate  deformation.  Furthermore,  the  shear,  twist,  and  rotary  inertia  effects  of  the  plate  are  taken  into 
consideration.  In  Eq.  1,  the  displacement  field  depends  linearly  on  the  coordinate  z  since  it  is 
assumed  that  the  plate  thickness  is  small  compared  to  the  other  two  dimensions. 
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3. 


NODAL  COORDINATES  AND  SHAPE  FUNCTIONS 


The  assumed  displacement  field  of  Eq.  1  includes  48  unknown  polynomial  coefficients,  which 
means  that  the  same  number  of  nodal  coordinates  is  required.  Figure  1  shows  the  element  used  in 
this  investigation  which  has  nodes  A,  B,  C  and  D.  In  Fig.  1,  a,  b,  and  t  are  the  length,  width  and 
thickness  of  the  plate  element,  respectively.  In  the  model  developed  in  this  investigation,  12 
coordinates  are  used  for  each  node.  The  coordinates  e7  of  a  node  j  can  be  chosen  as  follows: 


e7  = 


_  T  .  T  „  T 

T  drj_  drj_  drj_ 
J  dx  dy  dz 


(2) 


or :  or : 

where  the  vector  r,-  defines  the  global  position  vector  of  node  j  and  the  three  vectors  — — —  and 

dx  dy 


dr j 

— —  define  the  global  displacement  gradients  of  node  j.  Nodal  coordinates  of  one  element  can  then 
dz 


be  given  by  the  vector  e: 


e  = 


T  T  T 

eB  eC  eD 


(3) 


The  polynomial  coefficients  (a„  b ,•  and  c,)  can  be  replaced  by  the  nodal  degrees  of  freedom  [1 1],  and 
the  global  position  vector  r  can  be  rewritten  as  follows: 


(4) 


where  S  is  the  matrix  of  the  global  shape  function.  The  matrix  S  can  be  written  as: 


(5) 


S  =  [s{l  S21  S3 1  S4 1  S5l  S6 1  S-,1  Ss I  S9l  Sl0 1 

5„I  sn I  5„ I  514I  515I  51#I] 

where  I  is  a  3  x  3  identity  matrix  and 

5,  =  -(£  -1)(>7- l)(V  -  n  +  2£2  - 1  - 1  .  j ,  52  =  - 1)2  (»7  - 1)  , 

53  =-^(7-l)2(^-l)  ,  54  =f^-l)(^-l),  S5  =  (7-l)£(V  -7-3£  +  2£2), 

56  =-^2(^-l)(/7-l),  57  =6^(^-l)2,  S8=-ttf(7 7-1), 

59  =  -  3£  -  3/7  +  2r/2  +  2£“j,  Sl0  =  a£2rj(g  - 1),  $„  -6?72^(^-l) 

5,2  =  /#7,  5n  =7(#-l)^2-^-3^  +  2^2j,  514  =a^(^-l)2 

5,5  =  -b n  V - 1)  in- - 1) ,  5,6  =  - 1) . 

In  the  preceding  equation,  £  =  x/ci ,  r\  =  y/b  and  C,  -zlt . 

4.  RIGID  BODY  MOTION 

If  the  plate  element  undergoes  an  arbitrary  rigid  body  motion,  the  vectors  of  nodal  coordinates  (Eqs. 
2  and  3)  can  be  expressed  as  follows  [13]: 

r  ir 

ca  =  s.x  Ry  R=  Ai  A\  A i  A2  A2  A2  A3  As  A3  > 

r  ir 

eB  =  R.x+Aia  Ry  +  Aia  R:  +  A\a  Ai  Ai  Ai  A2  Ai  Ai  As  As  As  > 
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(6) 


-c 


Rx  +  | Cl  +  Aj 2b  Ry  "t"  A2\Q  ■*"  A22b  "t-  -^3 j#  +  A32b 


nr 


^11  ^21  ^31  ^12  ^22  ^32  ^13  ^23  ^33 


eD  ~ 


/?T  +  Ai2b  Ry+A22b  R-  +  A32b  An  A21  A3l  Au  A22  A32  Al3  A23  A33 


where  Atj  is  the  ij  th  element  of  the  3  x  3  rotation  matrix,  Rx  Ry  and  Rz  are  components  of  the 
position  vector  that  defines  the  global  translation  of  node  A  of  the  plate.  If  rotations  are  defined  in 
terms  of  three  Euler  angles  </) ,  9  and  y/ ,  the  rotation  matrix  takes  the  form: 


A  = 


cos  y/ cos (j)  -  cos  9  sin  <j>  sin  y/ 
cos  y/  sin  (j)  +  cos  9  cos  (j> sin  y/ 
sin  9  sin  y/ 


-  sin  y/cos(j>-  cos  9  sin  <j)  cos  y/ 

-  sin  y/  sin  <j)  +  cos  9 cos  (j)  cos  y/ 

sin  9  cos  y/ 


sin  9  sin  (j) 
-sintfcos^ 
cos  9 


(7) 


Substituting  the  expressions  of  nodal  coordinates  (Eq.  6)  into  Eq.  4  that  defines  the  global  position 
vector  at  an  arbitrary  point  gives: 


r  =  Se  = 


Rx  +  x(cos  y/ cos  <f>  -  cos#sin^sin^)-y(sin^cos^  +  cos0sin^cos^)  +  rsin#sin^ 

Ry  +  x(cos  i/ sin  <j>  +  cos  9 cos  ^sin  y/ )  +  y(cos  9  cos  (/)  cos  y/  -  sin  iy  sin  <p)  -  z sin  9  cos  (j)  , 
R.  +  .vsin#sin^  +  ysin0cosy  h-zcos# 


(8) 


which  demonstrates  that  an  arbitrary  rigid  body  motion  can  be  described  using  the  shape  function 
introduced  in  Eq.  5.  It  is  clear  from  the  analysis  presented  in  this  section  that  the  spatial  coordinates 
x  andy  represent  orthogonal  parameterization  of  the  mid  plane  of  the  plate  in  the  case  of  rigid  body 

dr  dr 

motion.  This  is  not  the  case  when  the  plate  deforms  since  in  this  case  —  and  —  do  not  remain 

dx  dy 


orthogonal. 
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5. 


PLATE  DEFORMATION 


By  using  global  slopes  instead  of  the  rotations,  no  assumptions  are  made  on  the  magnitude  of  the 
rotations  within  the  element.  This  allows  complex  shapes  to  be  represented  using  a  small  number  of 
elements.  As  an  example,  Fig.  2  shows  deformed  shapes  of  the  plate  when  one  element  based  on  the 
shape  function  of  Eq.  5  is  used.  The  deformation  in  Fig.  2a  is  obtained  by  using  the  following 
numerical  values  of  the  nodal  coordinates: 

e  =  [0.3  0.3  -0.2  1  0  1  0  0.5  1  1  1  1  1  0  0  1  0  1  0  1  -1  0  0  1 
0.7  0.7  -0.2  1  0  -1  0  1  -1  0  0  1  0  1  0  1  0  -1  0  1  1  1  1  l]r  . 

The  deformation  in  Fig.  2b  is  obtained  by  using  the  following  numerical  values  of  the  nodal 
coordinates: 

e  =  [0.3  0  0  1  0  1  0  1  0  1  1  1  0.7  0  0  1  0  1  0  1  0  0  1  1 

11010101000101010101011  i]r  . 

The  results  presented  in  Fig.  2  demonstrate  that  very  large  deformation  and  rotation  within  the  plate 
can  be  obtained  using  one  element.  This  important  isoparametric  property  is  very  crucial  in 
developing  efficient  shell  elements. 


In  the  case  of  rigid  body  motion  — ,  —  and  —  are  three  orthogonal  unit  vectors.  This  is 

dx  dy  dz 

no  longer  the  case  when  the  plates  deforms.  Figure  3  shows  the  norms  of  these  three  vectors  at 
points  of  coordinates  =  (£,  0.5,  0)  for  the  deformed  plate  shown  in  Fig.  2b.  The  dot 


are  shown  in  the  Fig.  4.  The  results  presented  in  this  figure 
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Sr  Sr  Sr 

show  that  the  three  vectors  — ,  —  and  —  are,  in  general,  no  longer  orthogonal  vectors  when  the 

dx  dy  dz 


plate  deforms. 


The  results  presented  in  Figs.  3  and  4  show  that  the  vectors  — ,  —  and  —  do  not  remain 

dx  dy  dz 

orthogonal  unit  vectors  as  the  plate  deforms.  The  coefficients  of  the  first  fundamental  form  of 
surfaces  can  be  used  to  shed  light  on  the  configuration  of  the  mid  surface  as  the  plate  deforms.  For 
the  mid  surface,  z  =  0  and  the  first  fundamental  form  in  this  case  is  defined  as  [7]: 

/  =  dr:= o  •  <7r._o  =  E  (dx)2  +  IFdxdy  +  G  (dy)2  ,  (9) 

where  E,  Fand  G  are  the  coefficients  of  the  first  fundamental  form  defined  as: 


E  =1  — 
dx 


dr 


F  = 

:=oydx):=0  ’  [dxj 


drY  f  dr^ 


.--o  [dy 


fdrY 


G  = 


J  r=0 


dyj^olfy):* 


dr 


(10) 


The  coefficients  E  and  G  measure  the  deviation  of  —  and  —  from  unit  vectors.  The  coefficient  F 

dx  dy 

Qr 

measures  the  deviation  of  —  and  —  from  orthogonal itv.  Figure  5  shows  these  three  coefficients 

dx  dy 


for  the  plate  shown  in  Fig.  2b  for  the  spatial  coordinates^,//,^)  =  (£,  0.5,  0). 


6.  PLATE  INERTIA 


In  the  absolute  nodal  coordinate  formulation,  the  mass  matrix  of  the  plate  element  can  be  obtained 
using  the  following  expression  for  the  kinetic  energy: 


(H) 
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where  p  is  the  mass  density  of  the  plate,  Fis  the  volume  and  r  is  the  absolute  velocity  vector.  The 
absolute  velocity  vector  is  linear  in  the  nodal  velocities  and  can  be  obtained  by  differentiating  Eq.  4 
with  respect  to  time  as  follows: 


r  =  Se  .  (12) 

Substituting  Eq.  12  into  the  expression  of  the  kinetic  energy  gives: 

T  =  ^erMe  ,  (13) 

where  M  is  the  mass  matrix  given  explicitly  in  the  appendix.  This  mass  matrix  is  constant  since  it 
only  depends  on  the  inertia  properties  and  dimensions  of  the  plate.  As  pointed  out  in  previous 
studies,  the  mass  matrix  obtained  using  the  absolute  nodal  coordinate  formulation  accounts  for  the 
effects  rotary  inertia  and  plate  twist  [9].  As  a  result  of  having  a  constant  mass  matrix,  the  centrifugal 
and  Coriolis  inertia  forces  are  equal  to  zero.  Having  a  constant  mass  matrix  allows  also  the 
development  of  an  efficient  algorithms  for  solving  the  nonlinear  equations  of  motion  of  multibody 
systems  [14,  18]. 

7.  ELASTIC  FORCES 

The  use  of  a  global  coordinate  system  for  the  definitions  of  the  nodal  coordinates  leads  to  a  simple 
expression  for  the  inertia  forces  and  a  non-linear  expression  for  the  elastic  forces.  Two  different 
methods  can  be  used  when  the  elastic  forces  are  derived  [13].  In  one  method,  a  local  element 
coordinate  system  that  is  used  for  the  description  of  the  element  deformation  is  introduced.  Another 
method,  which  is  used  in  this  study,  is  based  on  a  continuum  mechanics  approach  that  does  not 
require  the  use  of  a  local  element  coordinate  system.  The  latter  approach  leads  to  significant 
simplification  in  the  vector  of  elastic  forces  as  is  previously  shown  [3],  Furthermore,  by  employing 
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a  continuum  mechanics  approach  all  elastic  nonlinearities  are  taken  into  consideration  since 
nonlinear  strain-displacement  relationships  must  be  used. 


The  elastic  forces  of  the  plate  element  can  be  derived  by  using  the  virtual  work  or  the 
following  expression  of  the  strain  energy: 


v  =  e.dV  ,  (14) 


where  E  is  the  matrix  of  elastic  coefficients,  and  s,„  is  the  strain  tensor.  The  strain  tensor  can  be 
obtained  by  employing  the  Cauchy-Green  formula  [4]  as  follows: 


s„,=i(jrJ-l)  ,  (15) 


where  I  is  a  3  x  3  identity  matrix,  and  J  is  the  matrix  of  the  displacement  gradients.  Using  Eq.  4  the 
displacement  gradients  can  be  written  in  the  following  form: 


<9S  |  3S  i  dS  I 

— Le  — Le  — Le 

ox  dy  dz 

as,  as,  as, 

— -e  — -e  — -e 

dx  dy  dz 

as3  as3  as3 

— -e  — -e  — -e 

dx  dy  dz 


SLx.e  S^e  S,.e 

S2xe  S2ve  S2.e  ,  (16) 

^3.re  S3ve  S3re^ 


where  S,  is  the  /th  row  of  the  shape  function  matrix.  The  elastic  forces  of  the  plate  Q*  can  be 
obtained  by  differentiating  the  strain  energy  once  with  respect  to  the  nodal  coordinates,  that  is 


It  can  be  shown  that  the  shape  function  (of  Eq.  5)  leads  to  zero  elastic  forces  under  an  arbitrary  rigid 
body  displacement.  The  vector  of  the  elastic  forces  can  also  be  expressed  as  a  product  of  a  nonlinear 
stiffness  matrix  K  and  the  vector  of  nodal  coordinates  e  as  follows: 


Q*  =Ke 


(18) 


It  can  be  shown  that  the  stiffness  matrix  of  the  plate  element  can  be  written  as  follows: 


k=k,+k2+k3+k4+k5+k6 


(19) 


where 


/lie  bje - 1)  +  (A  +  IfA J b2e - ll  +  A ( efb3e - 1 


+ 


T  T 


e  b2e  - 1  (Ab,  +  (A  +  2//)b2  +  Ab 3)dV  , 


K2  =  ^  jb3  'lj^7b2e-lj  +  (A  +  2//)fe7'b3e-lj  +  A^e7'b1e-l 
erb[e-lj  (Ab2+(A  +  2//)b3  +  Ab,)t/F  , 


+ 


1  r  T 
K3  =-  lb, 

8,-!  1 

V 


T  \  f  T  \  f  T 

A [e  b2e-l  +  (/l  +  2 ju)  e  bje-1  +A  e  b3e-l 


+ 


c  b|  e  —  I  (/lb2  2/i)bj  4-  /lb2 ^jdV  , 


\  T  T  T  T  \  ~r  T  T  T 

K4=-//Jb4e  b4e  +  e  b4e b4  dV ,  K5  =  -/r  Jb5  e  b5e+e  b5eb5c/F, 

4  v 


|  y~  y*  y>  y 

K6=-iufb6-e  b6e+e  b6e  b  ^dV , 

A  J 
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where  /u  and  X  are  Lame’s  constant,  and 


bl  “  S3*Slz  +  S3^S2z  +  S3zS3z  »  b2  =  SUSLt  +  Sl^S2-v  +  - 

b3  =S2xS1>,  +82^82^  +S2zS3>),  b4  =SjTS1^,  +  S]yS2j,  +SlzS3j,  , 

b5  ~^l.T^lz +^l>>^2r +^lz^3z>  b6  =  ^2.t^1z  +  ^2>-^2z  +  ^2z^3z  • 

8.  IMPLEMENTATION  AND  NUMERICAL  RESULTS 

The  use  of  the  new  plate  element  developed  in  this  study  is  demonstrated  in  this  section  by  studying 
the  behavior  of  two  simple  pendulums.  Using  the  developments  presented  in  the  preceding  sections, 
the  equations  of  motion  of  the  finite  element  can  be  written  in  a  matrix  form  as  follows  [13]: 

Me  +  Ke  =  Qc.  ,  (20) 

where  Qe  is  the  vector  of  the  generalized  external  nodal  forces.  Using  the  definition  of  the  elastic 
forces,  the  preceding  equation  can  be  written  as: 

Me  =  Q  ,  (21) 

where  the  force  vector  Q  is  given  by 

Q  =  Qt  -Q*  (22) 

Since  the  mass  matrix  is  constant  the  equation  of  motion  can  be  efficiently  solved  for  the 
accelerations  using  the  following  equation: 
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e  =  M‘Q 


(23) 


The  pendulum  plate  models  in  this  investigation  are  shown  in  Fig.  6.  One  plate  model  has 
one  element,  while  the  other  model  has  three  elements.  In  the  two  models,  the  plate  is  connected  to 
the  ground  using  a  spherical  joint.  Therefore,  each  model  has  three  rigid  body  rotational  degrees  of 
freedom.  The  two  pendulums  are  assumed  to  vibrate  under  the  effect  of  gravity.  The  length,  width 
and  thickness  of  the  first  plate  are  0.3  m,  0.3  m,  and  0.01  m,  respectively.  The  material  of  the  plate 
is  assumed  to  be  isotropic  and  its  Young’s  modulus  is  1.0- 105  N/m2,  Poisson’s  ration  is  0.3  and 
mass  density  is  7810  kg/nv\  Figure  7  shows  the  motion  simulation  of  the  first  pendulum  and  Fig.  8 
shows  the  global  coordinates  of  point  H  (Fig.  6)  versus  time.  As  shown  in  the  Fig.  7,  one  element  is 
capable  of  representing  large  deformation  and  rotation  of  the  flexible  plate  under  the  effect  of 
gravity. 

The  second  pendulum  model  has  a  rectangular  plate  that  has  length,  width  and  thickness  of 
0.3  m,  0.6  m,  and  0.01  m,  respectively.  The  material  properties  of  the  plate  in  this  model  are 
assumed  to  be  the  same  as  the  properties  used  in  the  first  model.  Figure  9  shows  the  simulation  of 
the  plate  motion  under  the  effect  of  gravity.  Figure  1 0  shows  the  global  coordinates  of  point  H  (Fig. 

5r  (3i* 

6)  of  the  pendulum  versus  time.  The  components  of  the  vectors  — ,  —  and  —  along  line  /  (Fig.  6) 

dx  dy  dz 

on  the  mid  surface  of  the  plate  at  time  1.1  s  are  shown  in  Figs.  11,  12  and  13.  The  results  presented 
in  these  figures  demonstrate  that  the  global  displacement  gradients  at  the  nodal  points  are 
continuous. 

9.  SUMMARY  AND  CONCLUSIONS 

In  this  investigation,  a  new  plate  element  developed  using  the  non-incremental  absolute  nodal 
coordinate  formulation  is  presented.  The  general  plate  model  developed  in  this  investigation  can 
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describe  rigid  body  motion,  finite  rotations  and  an  arbitrary  large  deformation.  Continuity  of  all 
displacement  gradients  at  the  element  nodal  points  are  ensured,  thereby,  ensuring  the  smoothness  of 
the  mid-surface  of  the  plate  structures  at  the  element  nodes.  While  the  plate  element  developed  in 
this  study  relaxes  some  of  the  assumptions  of  the  conventional  and  Mindlin  plate  theories,  the  mass 
matrix  remains  constant  when  the  plate  undergoes  finite  rotations  and  experiences  large 
deformation.  It  is  demonstrated  in  this  investigation  that  complex  deformation  shapes  can  be 
modeled  using  small  number  of  elements.  By  using  global  slopes  instead  of  rotations,  no 
assumptions  are  made  with  regard  to  the  magnitude  of  the  rotations  or  the  deformation  within  the 
element.  A  continuum  mechanics  approach  is  used  to  obtain  the  plate  elastic  forces.  Nonlinear 
strain  -displacement  relationship  are  used  in  formulating  the  elastic  forces  and  as  a  result  the 
formulation  presented  in  this  paper  accounts  for  all  geometric  nonlinearities.  The  plate  element 
developed  in  this  study  is  of  the  isoparametric  type  and  can  be  used  to  develop  efficient  shell 
formulations. 
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APPENDIX 


The  mass  matrix  M  in  Eq.  13  is  defined  as  follows: 


Maa  Ma 

symm. 


The  submatrises  Maa,  and  can  be  written  as: 


- 


— —bml  JULpbO.l  —pbQ  I 
25200  12600^  “*  90^  • 

~^—b2ml  ~^—pb2Q\  —pb20\ 

315  200^  ^  60^ 

—pbl.\  — pbll 
105^  -  30^  * 

\pK' 


613  T  137  ,  ¥ 

- ml - bml 

12600  12600 

137  ,  T  1  ,2  T 

■—'“bml  -  ——b  ml 

12600  420 

Ho^1  -±f*QJL 

^pbQ} 

Mil  r  461  ,  T 

- ml - bml 

12600  25200 

~l—b2m\ 

315 


—pbOJ  —pbOl 
00y  ~z  360y  ~y 

jp4!e,i 

-ob2IJ  '-pbl,  I 


-^-pbO.l  —pbO  I 

2600  ”  —  90H  - 

~—pb2Ol  ~—pb2I  \ 

200y  ~  60^  - 

—pbl.. I  —pbl  l 

105^  “  30^  * 

l<*‘„  i 


M  of 


— — bml 
6300 

— pbO.l 
3150^  — 

hpb&1 

613  , 

- ml 

12600 
199  .  T 

-25200*"" 

i^1 


- pb2Q\ 

450^ 

-±-pb2Qv\ 
180^  - 

199  u  i 

- bml 

25200 

— 62otI 
630 

-355^ 


~jkpm 

-450^1 

— L  pbi  ,i 
280^  - 

-izpblj. 
90^  - 
137 

6300  Pb~: 

^r:pb20.l 

300  ' 

— -pb2I..\ 
\A0y  - 


ISO^1 


1)^1 

_  —pb20  i 
120^  ' 


613  , 

- ml 

12600 

199  t,  r 

- bml 

25200 

197  w 

- ml 

12600 

29  l  r 

- bml 

6300 

29 

—pbO  I 
3150^ 

— p&OvI 
45^  ~v 


25200 

— —  b2ml 
630 

So*2*1 

ilo^1 

29  r  I 

- bml 

6300 


-4t/?62O.I 

450^ 

— pb2IJ 

I  on"  >- 


^/o6ai  TTtpbOJ 
KT  360  ~ 


— L  pb2I_\ 

140^  “ 

-hpbIJ 

— —pbQ.l 

3150^  ** 

-Lpb2Q\ 


TopbI^ 


T5phI«l 

>4i 

■7^pt>% 1 


16 


M bb  = 
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where  I  is  a  3x3  identity  matrix,  /w  is  the  element  mass  and  p  is  the  mass  density  of  the  element  and 


~  12^,  1  ? 
0v=-at  ,  Q.  =~Q  t  ,  iyy 


1  3  r  1  3,  r  1  2,2 

-ot  ,  1~  =-a  t,  =-a  t  . 
j  j  4 
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Fig.  3.  Norms  of  vectors  — ,  —  and  — . 
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Fig.  4.  Dot  products 
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Fig.  6.  Pendulum  models. 
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Fig.  7.  Motion  simulation  of  the  first  plate  model. 
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Fig.  8.  Global  coordinates  of  point  Hof  the  first  pendulum 


t=0.1  s 


Fig.  9.  Configurations  of  the  second  plate  model. 
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Fig.  10.  Global  coordinates  of  the  free  end  (point  H)  of  the  second  pendulum. 
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Displacement  gradients  Displac 
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Fig.  11.  Components  of  the  vector 
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Fig.  12.  Components  of  the  vector 


